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A nonrandom complex method is described for application to optimization problems characterized 
by nonlinear objective and constraint functions involving continuous and/or discrete optimization 
variables. The method is a mutation of the "complex" method (involving a pseudo-random process) 
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1 . Introduction 

There are many engineering problems which have 
an unlimited number of adequate solutions, but which 
have an optimum solution with regard to some cri- 
terion such as weight or cost. This paper describes a 
nonrandom complex method for seeking the optimum 
solution for a class of such problems characterized by 
nonlinear objective and constraint functions of con- 
tinuous and/or discrete variables. The nonrandom 
complex method is a mutation of the 'complex" 
method of Box [l]. 2 Both the nonrandom method and 
the Box method have an important advantage over 
many other methods for nonlinear optimization in that 
they do not require computation of partial derivatives 
of the functions. 

The nonrandom method and the Box method begin 
a search by generating a configuration of points (an 
initial complex) whose coordinates are values of the 
independent optimization variables. At least one point 
of the initial complex (the initial point) must be known, 
in advance, to fall inside all constraint bounds of the 
problem. All other points of the complex are computed 
to fall inside or on the constraint bounds. 

In the Box method, the points of an initial complex 
are generated by a pseudo-random process. Numerical 
results by Hilleary [2], using the Box method, demon- 
strate that the rate of convergence to a solution is 
dependent on the character of the initial complex, and 
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2 Figures in brackets indicate the literature references at the end of this paper. 



that a poor initial complex may result in failure of the 
method to converge. 

In the nonrandom complex method reported here, 
the character of the initial complex is influenced by 
the initial point to a much greater degree than in the 
Box method. Thus, a good initial point leads to a rela- 
tively good initial complex. This characteristic is used 
to advantage by generating a new initial complex about 
a new initial point if the convergence of a search se- 
quence becomes excessively slow or stops. The same 
characteristic also facilitates a more efficient discrete 
variable optimization than would be expected from the 
Box method. The nonrandom method imposes a lower 
bound on initial complex size by selecting initial com- 
plex points in such a way that at least one-half the 
constraint bounded range of each independent variable 
(in the one-dimensional subspace containing the initial 
point) is included within the complex span. 

The nonrandom complex method has been used 
successfully for the following two minimum weight 
structural analysis problems: 

1. An ^-degree elliptical elastic ring with sinusoidally 
varying cross-section dimensions, loaded by two forces 
acting in opposite directions along a diameter; and 

2. A rib-stiffened simply supported plate, loaded by* 
in-plane compressive forces. 

The ring problem has six continuous independent 
variables, and the plate problem has four independent 
variables, two of which are discrete. Both problems 
have nonlinear objective and constraint functions. A 
large class of optimization problems fall within the 
range of mathematical complexity represented by 
these two problems. 
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List of Symbols 

%k independent variable of the optimization 

problem, 

gi implicit constraint function, 

hi implicit constraint lower bound, 

ki implicit constraint upper bound, 

E m m-dimensional Euclidean space, 

Pj a point of the optimization complex, 

ah lower explicit bound on Xk, 

bk upper explicit bound on Xk, 

a reflection factor of complex optimization method, 

m the number of independent variables of the 

optimization problem, 

p the number of independent discrete variables 

of the optimization problem, 

L total force applied to ring, 

P =L/2, force applied to one quadrant of ring, 

8 curved bar deflection, one-half ring deflection, 

r, 6 polar coordinates of ellipse, 

n degree of ellipse, 

a, b semimajor and semiminor axes of ellipse, 

h ring cross-section width, 

t ring cross-section thickness normal to tangent 

of ellipse, 

c average of maximum and minimum ring width, 

d average of maximum and minimum ring thick- 
ness, 

e dimensionless parameter of ring width variation, 

/ dimensionless parameter of ring thickness 

variation, 

s arc length variable measured along ellipse, 

M moment, 

N normal force, 

V shearing force, 

k curvature of ellipse, 

g integral defined in equations (7) and (8), 

v Poisson's ratio, 

E modulus of elasticity, 

a normal stress, 

S design maximum stress, 

A scale factor, 

W weight, 

y material density, 

Ni normal force per unit width of plate, 

/, w length and width of rib-stiffened plate, 

t p flat plate thickness, 

t r rib thickness, 

b p rib spacing, 

b r rib width, 

H =t p -\-b r , total thickness of rib-stiffened plate. 



2. The Nonrandom Complex Method 



The following paragraphs describe the nonrandom 
complex method as it is applied to a continuous vari- 
able problem. Adaptation of the method for discrete 
variable optimization is described later. 

The nonrandom complex method searches for a set 
of values of the independent continuous variables 



#i, #2 9 • • ., %m which gives a minimum value of an 
arbitrary objective function f(xi, x%, • • .? #m), subject 
to explicit constraints on the variables, and subject 
to implicit constraint bounds on arbitrary functions of 
the variables. The implicit constraint functions 
gi(xi, X2, . . ., x m ) are subject to inequalities of the 
form hi ^ gi ^ ki, where the bounds hi and ki represent 
arbitrary limits and either one or both may appear. 

To begin the method it is required that a point 
{x[, #2, . . ., x' m ) , in m-dimensional Euclidean space, 

E m \ be known which is feasible (i.e., does not violate 
any constraint) and which does not lie on any constraint 
bound. Using this initial point, a set of 2m +1 points 
in E m is then generated by the following procedure: 

1. Denote the jth point in the set of 2m+ 1 points by 
1 j . j_j e 1 1 1 1 x ■* , x 9 9 • • * 1 x ni j • 

2. Denote the explicit constraints on the independent 
variables by a* ^ Xk ^ bk for k= 1 , 2, . . . , m. Set k = 1. 



3. Let P2A— (*i, x{ 



2-> 



' ) Xk-\, Ok-> Xk 



+ 1? 



Xm ) 



This point is the same as the initial point P\ except in 
the A:th coordinate, where x' k has been replaced by 
its explicit upper bound bk- 

4. Evaluate all implicit constraint functions at Pzk- 
If all are satisfied, move on to step 5. If any implicit 
constraint is violated, set 



Pn = 



b* + 



X 1 



•\? 1 - -A> 



2' 



X 



A'-l' 



X 



A:+l> 



• * -A. 



m 



This moves P>ik one-half the distance toward the initial 
point P\. Once again, evaluate the implicit constraint 
functions at P>k, go to step 5 if they are satisfied, or 
move P 2 k one-half the remaining distance to P\ other- 
wise. Repeat the process until P>zk is feasible. 

O. l_jet 1 2A*+ 1 " " \X 1 , X>2 1 • • • 1 Xj._^ , (Ik 1 -^A*+l ' ' * ' ' ni )' 

This point is the same as the initial point P\ except in 
the A:th coordinate, where x' k has been replaced by 

its explicit lower bound dk. 

6. Proceed with point P2A+1 as in step 4, moving 
P2A+1 one-half the distance to the initial point Pi 
until P2A+1 is feasible. 

7. Reiterate steps 3 through 6 for A: = 2, 3, . . . , m. 
This gives a collection of 2m +1 points {Pi,P 2 , • • ., 
Pzm+i) which is the initial complex. 

The method proceeds to search for an optimum by 
repeatedly altering the complex, one point at a time, 
as follows. The objective function f(x\ , X2, . . ., x m ) 
is evaluated at each complex point Pj. The complex 
point with the worst (greatest) objective function value 
is selected, and the centroid of the remaining 2ra points 
is computed. The worst complex point is then discarded 
and replaced by a new point, with the exception that 
the second worst point is discarded if the worst point 
is also the newest complex point. That is, a given com- 
plex point is not permitted to be changed in consecu- 
tive iteration cycles. A provisional replacement point 
is computed which is: 

1. Collinear with the discarded (worst or second 
worst) point and the centroid of the 2m complex points 
excluding the worst point; 
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2. On the side of the centroid opposite the discarded 
point; and 

3. Located a times as far from the centroid as was 
the discarded point, where a is a constant parameter 
called the reflection factor. 

If the provisional replacement point violates an explicit 
constraint on any variable, the corresponding coordi- 
nate is set equal to the violated constraint bound. If 
the new provisional point then violates an implicit 
constraint, it is repeatedly replaced by a new point, 
located one-half the distance to the centroid of the 2m 
points, until a feasible new complex point is found. 
The new centroid of all 2m + 1 complex points is then 
determined, and the objective function at that cen- 
troidal point is evaluated and used as an index of con- 
vergence. This completes one cycle of the iterative 
search process. 

During a search, a series of cycles may yield an 
oscillating value of the centroidal objective function 
with little or no net convergence. If a preassigned 
number of consecutive oscillating cycles yield no net 
improvement in the centroidal objective function value, 
a new initial complex is generated, using the best 
centroid yet computed as the initial point. 

A search is stopped if a preassigned number of 
consecutive cycles yield no change in the value of the 
centroidal objective function greater than a preas- 
signed small amount. This condition may indicate 
that an approximation of the global optimum has been 
found. But, the same condition could indicate that a 
relative (nonglobal) optimum has been found, or that 
the complex has collapsed to a small size due to the 
reflection factor being too small or due to some other 
characteristic of the particular problem. Cycle count 
and elapsed time stops are also used to prevent exces- 
sive computer time. 

A search may also be stopped if the centroid of the 
2m nonworst complex points lies on or outside a con- 
straint bound. This can be caused by a nonconvex 
constraint surface or, if the complex is small, it can 
be caused by numerical rounding errors. 

If a search is stopped before a satisfactory approxi- 
mate optimum is found, it can be restarted by generat- 
ing a new initial complex, using one of the better 
complex centroids already computed as the initial 
point. A centroid somewhat offset inside all constraint 
bounds, and not necessarily the centroid with the 
lowest objective function value, is used for this 
purpose. 



3. Discrete Variable Optimization 

If variables Xi through x n of the m optimization 
variables, where p ^ m, are limited to known sets of 
discrete values, the following adaptation of the non- 
random complex procedure is applicable. All m 
optimization variables are temporarily assumed to be 
continuous and an optimum continuous variable 
solution x— (%u %2, • . ., x m ) is found by the non- 
random complex method. A set of 4? neighboring dis- 



crete value points is then generated by the following 
procedure: 

1. Select the two discrete values of the A:th discrete 
variable, where k ^ p, that are greater than Xk and 
that most nearly approximate xk* 

2. Select the two discrete values of x^ that are less 
than or equal to Xk and that most nearly approximate 

Xk- 

3. Reiterate steps 1 and 2 for A:=l, 2, . . . , p. This 
gives 4 discrete approximate values for each #*, 
k^p. 

4. Generate the 4> p distinct combinations of the dis- 
crete values obtained by step 3. Each combination is a 
set of the discrete coordinates for a discrete approxi- 
mation to the continuous variable solution point 
x = (xu x~2, . . ., x m )- The remaining coordinates for 
the 4^ discrete points that surround and approximate 
x are the continuous variable coordinates Xk, where 
k > p. 

The objective and constraint functions are then 
evaluated at each of the 4 ;> discrete value points, and 
the best (lowest objective function value) and second 
best feasible points are selected. These two points 
are used as the initial points for two nonrandom 
complex optimization searches. For each search, the 
discrete variables are held constant and the continuous 
variables are optimized within the subspace defined 
by the discrete coordinates. The two resulting solutions 
are then compared. If the better of the two solutions 
results from the better initial point, that solution is 
assumed to be the optimum solution to the problem. 
If the better solution results from the second best 
initial point, the next (third) best feasible point of the 
4^ discrete points is used as the initial point for a third 
continuous variable optimization within the discrete 
coordinate subspace containing the point. The process 
is repeated similarly until a solution results that is 
worse than a previous solution. 

If the number of independent discrete variables is 
such that 4^ discrete value points would be an exces- 
sively large number to evaluate, an alternative pro- 
cedure is advisable. A reasonable compromise would 
be to generate and evaluate the 2 P adjacent discrete 
value points that bound the continuous variable 
solution x. 



4. Elastic Ring Problem 

An earlier paper [3] reported a minimum weight 
structural analysis of an n-degree elliptical elastic 
ring, with sinusoidally varying cross section dimen- 
sions, by a sequential grid method. The same elastic 
ring optimization problem is here used to demonstrate 
application of the nonrandom complex method to a 
difficult continuous variable problem. 

The objective of the optimization problem is to find 
the dimensions of the minimum weight ring that would 
satisfy the arbitrary force capacity, deflection, and 
dimensional requirements of the load supporting ele- 
ment of a force transducer. The force capacity require- 
ment is satisfied by a scaling procedure. King deflection 
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2P = APPLIED FORCE 




28 = CHANGE IN 
THIS DIMENSION 




K 



!Lfl4aibL_ 



2P 



Figure 1. 



A 1,000,000 lb capacity n-degree elliptical ring with 
sinusoidally varying thickness. 



and critical ring dimensions are expressed as implicit 
constraint functions. 

An elastic ring force transducer is loaded by two 
forces acting in opposite directions along a diameter, 
as shown in figure 1. Corresponding ring deflections 
28 are measured and related to applied forces 2P 
by a calibration factor. A structural analysis (and 
FORTRAN computer program) for the elastic ring was 
given in the earlier paper [3]. The ring was analyzed 
as a thick curved bar, assuming small deflections, 
plane strain, Winkler-Bach bending stress distribution, 
uniform normal stress distribution, and parabolic 
shear stress distribution. The following paragraphs 
summarize the structural analysis and give the equa- 
tions required for the nonrandom complex optimization. 

For a ring of the shape shown in figures 1 and 2, 
the locus of centroids of ring cross sections is an 
M-degree ellipse defined by the equation 



r=ab(b n cos n 0-\-a n sin" 6Y lln , n^2. 



(1) 



Ring rectangular cross-section width h and thickness 
t are given by 



h = c(l — e cos 20), 



(2) 



and 



; = d(l-/cos2<9), 



(3) 



in which c and d are positive parameters and e and/ 
are parameters of absolute value less than unity. 



Mjr/2 



Vt/2 




FIGURE 2. Loads acting on one quadrant of a ring. 



Using these equations, the shape of a ring within this 
geometric class can be specified by the seven param- 
eters a, b, c, tZ, e,/, and n. If a = b, n = 2, and e=f= 0, 
the ring shape is circular with uniform rectangular 
cross section. 

Figure 2 shows the loads acting on one quadrant of 
the ring. The resultant moment and forces acting on a 
typical normal cross section are 



M = M v ,2 ~ Pr cos 6 = PM' 



(4) 



N=P ( sin 6 ^+r cos 6 jf) = ^'> (5) 



and 



V=P ( r sin 6 ~ cos 6 ^)=PV, 

as as 



(6) 



in which rfs is the increment of arc length of the n- 
degree ellipse. With the assumed stress distributions, 
the strain energy of the ring quadrant is obtained by 
integrating over the quadrant volume. Castigliano's 
strain energy theorem is used to determine the mo- 
ment and deflection at the point of load application 
to be 



M wt2 = P 



I 2 {1 -v 2 )gr cos , f^ 2 k sin 6 . 

as — — : dr 



h 



ht 



n l 2 kr cos 
ht 



dO 




' 2 (1 - v 2 )g 
h 



ds 



(7) 



and 



. _p r-/2 i 

K,2 ~E Jo ht 



/V' 2 + ^ (1 + v)V' 2 - (1 - v*)gtr cos 6M' + kM'N' - kr cos ON' 



ds, 



(8) 
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in w 



hich 



k 



^■S)'-'"" 



Pr 



2 



r* + 



de, 



3/2 



£ 



J/ 



< 




t 


Y 




(l + ky) 


I 


y* 


dy 




(l + ky) 



1 2 



dy, 



^=Poisson's ratio, 

E= modulus of elasticity, 

and y is the distance from a cross-section centroid to 

a point in the cross section, positive in the outward 

direction. The stress in the ds direction on the inner 

surface of the ring is given by 



cr T 



M'k 

ht 



M't 



h(2-kt) 



f 



+ 



\ (l + ky) 



dy 



K 

ht 



(9) 



For the optimization analysis, it is assumed that the 
force resisting capacity of the ring is limited by the 
maximum value of oy for all 6. The possibility of capac- 
ity being limited by some other component of stress 
or combined stress should be checked in any final 
design analysis. 

To make each numerical solution have the required 
force capacity, all length dimensions in the r-6 plane 
are scaled by a factor X that will make the maximum 
value of ctt, for the entire ring, equal a prescribed 
maximum design stress S. Cross-section width h is 
not scaled. The scale factor is 



X 



2S 



M'k 
ht 



M't 



h(2-kt) 



l 



f 



+ 



(l + ky) 



dy 



ht 



, (10) 



max. 



in which S = design maximum stress, 
L\ = required force capacity, 
and the maximum value of the bracketed function for 
all 6 is used. 

In the optimization analysis the above equations are 
first solved for a ring of unit value of mean radius, 
that is 

^=i. (id 

Since a and b are functionally dependent, the ring 
shape can be specified by the six independent shape 
parameters a, c, d, e, /, and n. After the structural 
equations are solved for the unit size ring (using 
Simpson numerical integration), the results are scaled 
to full size using the scale factor X. The width param- 
eter c is not scaled during the entire analysis, and the 
shape parameters e, /, and n are nondimensional and 
are not scaled. 

The six shape parameters a, c, d, e,/, and n are the 
six independent variables xi, X2, . . . , x$, respectively, 
of the constrained optimization problem. The objec- 



tive function /(^i, x 2 , . . ., x 6 ) is the scaled ring weight 



JF = 4y\ 2 I"' ht ds, 
o 



(12) 



in which y is the material density. The following six 
implicit constraint functions g\ are applied: 



#i = 2X8 tt/2 = total ring deflection, 



(13) 



g 2 = Xd(l + |/| ) = maximum thickness of ring 

cross section, (14) 

g 3 =\(2a-\-d+ df) = outside width of ring, (15) 

g 4 = X(4 — 2a + d-\- df) = outside height of 

ring, (16) 

g$ = c(l + \e\) = maximum width of ring cross 

section, (17) 

and #6 = X (4 — 2a — d— df) = inside clear height of 

ring. (18) 

Implicit lower bounds hi and/or upper bounds /c, are 
imposed on the functions #* to limit ring deflection and 
dimensions. 

The scaling procedure used here could have been 
avoided by expressing force capacity as an implicit 
constraint function with a lower bound. It is not known 
whether such a procedure would be more efficient than 
the scaling procedure used. 

Numerical solutions were obtained for seven dif- 
ferent minimum weight elastic ring problems. Tables 
1 and 2 give the initial points, the constraint bounds 
and the solutions obtained. The material properties 
were: S= 150,000 lb/in 2 , E = 30,000,000 lb/in 2 , v=0.3, 
and 7=0.29 lb/in 3 . The computations were coded in 
FORTRAN V and run on a UNIVAC 1108 computer. 
Time required for a solution ranged from 4 min to 
18 min with a mean of 9 min. The number of initial 
complex regenerations required for a solution ranged 
from 2 to 21 with a mean of 9. A reflection factor of 
1.5 was used for all seven solutions. 

For each of the two 1,000,000 lb capacity ra-degree 
elliptical rings, the best solution reported in the earlier 
paper [3] (weight = 539.08 lb) was used as the initial 
point. For the other five problems, the initial feasible 
points were obtained by using the relatively simple 
equations for a thin circular ring of uniform cross 
section. The thin circular ring equations, although they 
may be in error by several percent, give a good enough 
approximation of ring characteristics for beginning an 
optimization search. The coordinates obtained by the 
thin ring equations were offset slightly inside their 
explicit bounds, where that was necessary, to comply 
with the requirements of the optimization method. 
To hold a particular parameter constant, the upper and 
lower explicit bounds were set equal to the initial 
point coordinate value. 
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Table 1. Initial points and constraint bounds for elastic ring problems 






Force 


Problem 
dimension 












Initia 


point coordinates 


and explicit bo 


unds 
















Implicit bounds 


capacity 


a\ 


a 


6, 


a 2 


c 


b 2 


a 3 


d 


63 


</ 1 


e 


64 


a 5 


/ 


h 


a 6 


n 


be 


h x 


k x 


h 


k 3 


&4 


h 


A. 


Lb 




Inch 


Inch 


Inch 


Inch ^ 


Inch 


■ 

Inch 


Inch 


Inch 


Inch 


. 
















' 


Inch 


Inch 


Inch 


Inch 


Inch 


■ 

Inch 


Inch 


100,000 


5 


1.0 


1.01 


1.3 


3.0 


5.0 


6.0 


0.05 


0.151 


0.3 


0.0 


0.0 


0.0 


-0.9 


0.0 


0.9 


2.0 


2.1 


20.0 


0.25 


0.55 


4.0 


40.0 


40.0 


6.0 


9.0 


300,000 


5 


1.0 


1.01 


1.3 


4.0 


8.0 


10.0 


.1 


.188 


.4 


0.0 


0.0 


0.0 


-0.9 


0.0 


0.9 


2.0 


2.1 


20.0 


.25 


.55 


4.0 


40.0 


40.0 


10.0 


9.0 


1,000,000 


2 


1.0 


1.00 


1.0 


6.0 


12.0 


14.0 


.1 


.250 


.5 


0.0 


0.0 


0.0 


0.0 


0.0 


0.0 


2.0 


2.0 


2.0 


.25 


.55 


7.0 


50.0 


50.0 


18.0 


9.0 


1,000,000 


5 


1.0 


1.21 


1.3 


6.0 


13.3 


14.0 


.1 


.305 


.5 


0.0 


0.0 


0.0 


-0.9 


0.5 


0.9 


2.0 


6.0 


20.0 


.25 


.55 


7.0 


50.0 


50.0 


18.0 


9.0 


1,000,000 


6 


1.0 


1.21 


1.3 


6.0 


13.3 


14.0 


.1 


.305 


.5 


0.9 


0.0 


0.9 


-0.9 


0.5 


0.9 


2.0 


6.0 


20.0 


.25 


.55 


7.0 


50.0 


50.0 


18.0 


9.0 


3,000,000 


5 


1.0 


1.01 


1.3 


14.0 


18.0 


20.0 


.25 


.311 


.5 


0.0 


0.0 


0.0 


-0.9 


0.0 


0.9 


2.0 


2.1 


20.0 


.25 


.55 


9.0 


60.0 


60.0 


20.0 


9.0 


9,000,000 


5 


1.0 


1.01 


1.3 


20.0 


28.0 


32.0 


.3 


.387 


.5 


0.0 


0.0 


0.0 


-0.9 


0.0 


0.9 


2.0 


2.1 


20.0 


.25 


.55 


15.0 


85.0 


85.0 


32.0 


9.0 


1 




































































Table 2. 


Solutions for elastic ring problems 














Force 


Problem 
dimension 




Ring s 


.hape parameters of unit size ring 






Implicit constraint function values 




Scale 
factor 


Ring 

weight 


Cycles 


capacity 






























a 


c 


d 


e 


f 


n 


g\ 


#2 


g3 


g4 


gh 


g6 








Lb 




Inch 


Inch 


Inch 








Inch 


Inch 


Inch 


Inch 


Inch 


Inch 




Lb 




100,000 


5 


1.07509 


5.35807 


0.14689 


0.00000 


0.79775 


5.36445 


0.2502 


1.5078 


12.4472 


12.0703 


5.3581 


9.0546 


5.710 


52.40 


210 


300,000 


5 


1.09246 


8.40706 


.19659 


.00000 


.67392 


4.52213 


.2503 


2.3088 


15.7788 


15.0431 


8.4071 


10.4256 


7.016 


163.17 


201 


1,000,000 


2 


1.00000 


12.32346 


.27146 


.00000 


.00000 


2.00000 


.2501 


4.1548 


34.7653 


34.7653 


12.3235 


26.4558 


15.305 


1427.90 


190 


1,000,000 


5 


1.19033 


11.87379 


.31110 


.00000 


.56132 


9.39859 


.2509 


3.8990 


20.2053 


16.8975 


11.8738 


9.0996 


8.027 


492.49 


179 


1,000,000 


6 


1.19717 


11.77591 


.31561 


- .02669 


.56049 


8.15270 


.2512 


3.9900 


20.5215 


16.9984 


12.0902 


9.0184 


8.102 


497.78 


263 


3,000,000 


5 


1.02028 


16.12696 


.33864 


.00000 


.62951 


5.93016 


.2504 


6.6805 


26.2232 


30.4025 


16.1270 


17.0414 


12.107 


1686.69 


522 


9,000,000 


5 


1.02117 


25.64302 


.45811 


.00000 


.51622 


5.69466 

. 


.2503 

he 


9.7050 


31.6321 


37.0573 


25.6430 


17.6474 


13.972 


4819.54 


414 



For the five n-degree elliptical rings with uniform 
cross section width, the ratio of weight to force capac- 
ity ranges from 4.92 X 10~ 4 to 5.62 X 10 4 . At least a 
rough proportionality between weight and force capac- 
ity might be expected for structures of this geometric 
class. Thin ring equations would predict a constant 
ratio of weight to force capacity for uniform cross 
section, circular rings having the same deflection at 
capacity load. 

In each of the seven solutions, ring deflection, gi, 
is near its lower bound, hi = 0.25 in. In three of the 
solutions, inside clear height, #6, is near its lower 
bound, hQ = 9.0 in. Exact minimum weight solutions 
would apparently lie exactly on at least one of the 
bounds. Continued iteration would have resulted in 
closer approximations to exact solutions. But the 
searches were stopped because of the relatively large 
computer costs involved and because the convergence 
rate decreases as an optimum is approached. 

The seven solutions include three 1,000,000 lb 
capacity rings of different geometric complexity. 
One ring is circular with uniform cross section; one 
is rc-degree elliptical with constant cross-section width 
and sinusoidally varying thickness; and one is ra-degree 
elliptical with sinusoidally varying width and thickness. 
The three cases are respectively two-dimensional, 
five-dimensional, and six-dimensional optimization 
problems. All three cases have the same implicit 
constraint bounds. Both n-degree elliptical rings have 
a weight of about 35 percent of the weight of the cir- 
cular ring, and the outside dimensions of the more 
complex rings are significantly less than those of the 
circular ring. Since the five-dimensional problem is a 
special case of the six-dimensional problem, an exact 
solution to the six-dimensional problem would have a 
weight equal to or less than any solution to the five- 
dimensional problem. Although both elliptical rings 
have about the same weight, the ring with nonuniform 
cross-section width would probably cost more to manu- 
facture than would the other elliptical ring. 



modes to be applicable to their problem. The follow- 
ing specialized forms of the buckling equations used 
in their study are assumed to apply to the present 
problem: 

1. Critical load for gross buckling of the simply 
supported, rib-stiffened plate, 



i V rr 
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(19) 



in which 



1 



3(1 -v 1 ) 




tj>V (It 



+ 



V 



2(1- -v 2 ) \ HI \H) \b 



i>> 



x 



1+^ 





D< = 




E= modulus of elasticity, 
v= Poisson's ratio, 



5. Rib-Stiffened Plate Problem 

The minimum weight analysis of a rib-stiffened 
elastic plate is presented here to demonstrate the 
application of the nonrandom complex method to a 
discrete variable problem. It is assumed that the rib- 
stiffened plate is built-up from flat plates that are 
available in only a limited number of thicknesses. 

Figure 3 gives a schematic drawing of a rib-stiffened 
plate of the type studied. The flat plate is reinforced 
by a rectangular grid system of ribs having the same 
dimensions and spacing in both directions. The rib- 
stiffened plate is simply supported along all four edges 
and loaded only with a compressive load Ni in the 
plane of the plate. 

Schmit, Kicher, and Morrow [4] have studied the 
application of gradient optimization methods to a 
more general class of rib-stiffened plates, subjected 
to combinations of both shear and biaxial normal 
loads. They considered three possible elastic buckling 



and the length dimensions are those indicated in 
figure 3. 

2. Critical load for local buckling of a stiffening rib 
as a flat plate, simply supported along three edges 
and free along the fourth edge, 



Ncr = 



^(irm[(^y+H- 

(20) 
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3. Critical load for local buckling of a square panel 
as a flat plate, spanning between two pairs of ribs and 
simply supported along all four edges, 
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(21) 



The critical material stress in the rib-stiffened plate 
is assumed to be equal to the applied load divided by 
the gross cross-sectional area of the flat plate and 
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Figure 3. Rib-stiffened plate. 



ribs. That is, 



CT/ = 



iVi 



tp + 6,- 



*r 



(22) 
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The gross weight of the rib-stiffened plate is given 
approximately by 



W=ylwt 



p 




(23) 



in which y is the material density. Equation (23) does 
not include the net weight resulting from the fastening 
of the stiffening ribs to the larger flat plate or to each 
other. 

For the particular numerical problem presented here, 
plate length and width, / and w, are held constant and 
the four plate dimensions t p , t r , b p , and b r are the vari- 
able parameters of the optimization search. The plate 
thicknesses, t p and t r , are discrete variables. The 
objective function is the gross weight of the plate, W. 
The following implicit constraint functions are applied: 



gi = N cr — gross buckling load, 
g-i — N cr — rib buckling load, 
g3 = N C r :== panel buckling load, 
g4 = (Ji — critical stress, 



(24) 
(25) 
(26) 

(27) 



and #5 = t p + b r = total stiffened plate thickness. (28) 

The applied load that the plate is required to support 
is a lower bound on g u go, and gs. The design stress, 
5, is an upper bound on^- 



A numerical solution was obtained for a rib-stiffened 
plate 30 in wide by 40 in long with an applied compres- 
sive load /V/=350 lb/in. Lower and upper bounds of 
0.5 in and 0.7 in, respectively, were imposed on total 
stiffened plate thickness, #5. Explicit constraint bounds 
on the variable optimization parameters were as 
follows (in inches): 



0.005 ^ t p ^ 0.500 



0.005 



2.000 



7> 



0.500 



6.000 



0.100^ b r ^ 0.700. 

The material properties were: 5 = 20,000 lb/in 2 , 
£"=10,500,000 lb/in 2 , *> = 0.32, ^and y = 0.101 lb/in 3 . 
Forty-four discrete plate thicknesses from 0.0040 in 
to 0.5800 in were assumed to be available. The compu- 
tations were coded in FORTRAN V and run on a 
UNIYAC 1108 computer. Total computer time re- 
quired for the solution was 55 s. A reflection factor of 
1.6 was used. 

It was determined by direct computation that a 
0.2043 in thick plate was the thinnest available plate 
that would support the applied load as a flat plate, 
without stiffening ribs. Such a plate would weigh about 
24.76 lb. The 0.2043 in flat plate thickness combined 
with a reasonable set of rib dimensions were used as 
the initial point coordinates for the first continuous 
variable search. 

Table 3 gives the initial points and results for the 
series of five optimization searches employed in finding 
a discrete variable solution. The result of the first 
continuous variable search was a solution point near 
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TABLE 3. Initial points and results for rib-stiffened plate problem 





Optimization 


Variable parameters 






Implicit constraint function values 




Gross 


Cycles 


search 


















weight 


























h 


tr 


b„ 


b r 


Ncr 


Ner 


N„ 


crt 


H 








Inch 


Inch 


Inch 


Inch 


Lb/Inch 


Lb/Inch 


Lbllnch 


Lb/Inch 2 


Inch 


IJb 






0.20430 


0.204™ 


4.00000 


0.30000 












28.3805 




First 


(Initial point) 














continuous 
























.04541 1 .04526 


3.44039 


.64625 


350.00 


945.27 


371.19 


6492.0 


0.69166 


7.5510 
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(Result) 






















.04767 


.05041 


3.72198 


.63748 


350.10 


1242.35 


365.24 


6216.3 


.6asi5 


7.R563 




Second 


(Initial point) 




















continuous 


























.03068 


.02% 1 


2.00111 


.64446 


350.10 


350.00 


374.79 


8703.1 


.67514 


6.0128 


215 




(Result) 






















.03068 


.02%1 


2.00110 


.64444 


350.10 


350.10 


374.91 


8702.2 


.67512 


6.0133 




Third 


(Initial point) 




















continuous 


* . 






















.02989 1 .02988 


2.000(H) 


.64409 


350.15 


350.51 


350.00 


8858.2 


.67398 


5.9375 


111 




(Result) 






















.03200 1 .03200 


2.00000 


.64409 


376.95 


430.72 


430.46 


8273.2 


.67609 


6.3565 




First 


(Initial point) 




















discrete 


1 






















.03200 | .03200 


2.19138 


.64465 


350.00 


406.67 


350.00 


8451.3 


.67665 


6.1436 


91 




(Result) 






















.03590 | .02850 


2.00000 


.64409 


351.38 


363.79 


575.24 


7764.3 


.67999 


6.56(H) 




Second 


(Initial point) 




















discrete 


1 






















.03590 1 .02850 


2.08008 


.65105 


350.01 


350. 00 


528.17 


7809.0 


.68695 


6.4985 


114 




(Result) 

1 





















the gross buckling constraint and not far from both the 
panel buckling constraint and a total plate thickness 
constraint. The continuous variable search was re- 
started, in a second computer run, using a complex 
centroid somewhat offset inside all constraints as the 
initial point. The result of the second continuous vari- 
able search was a solution point near a rib spacing 
constraint, near the gross buckling constraint, near 
the rib buckling constraint, and not far from the panel 
buckling constraint. The continuous variable search 
was again restarted, in a third computer run, and the 
resulting solution point was near a rib spacing con- 
straint and all three buckling constraints. As a con- 
tinuation of the third computer run, the two best 
discrete solution points near the continuous variable 
solution were selected and used as initial points for 
two searches within the subspaces defined by the 
discrete coordinates. The better of the two discrete 
solutions resulted from the better discrete initial point. 
Therefore, no additional discrete subspace searches 
were executed. The gross weight for the discrete solu- 
tion is about 25 percent of the weight of the thinnest 
available flat plate that would support the applied load 
without stiffening ribs. Initial complex regeneration, 
to overcome excessively slow convergence, was not 
required during any of the five optimization searches. 



6. Discussion 

Unfortunately, a general method is not available 
for determining whether or not a solution found by 
the nonrandom complex method is (approximately) 
the true optimum. But a careful study of the search 



record can yield a basis for judgment. The study 
should include not only the final values (centroidal) 
of the variables and functions, but also the varia- 
tion of these values for the entire search. The entire 
complex, at regular cycle intervals, should also be 
examined. 

Multiple searches along different paths can be used 
to test apparent optimum solutions. The search path 
depends on the objective and constraint functions, 
the constraint bounds, the initial point, the reflection 
factor, the complex regeneration criterion, and nu- 
merical errors. A change in one or more of these 
factors can cause a search to progress along a differ- 
ent path. Of course the rate of convergence is, in 
general, different for different search paths. 

Discrete variables can occur in practical problems 
in several different forms. The rib-stiffened plate 
problem involves two discrete variables with tabu- 
lar noninteger values. That problem could have been 
written to include an integer optimization variable by 
expressing rib spacing as a function of the number 
of equally spaced longitudinal ribs. Schmit and 
Mallett [5] have discussed the optimization of ma- 
terial selection by initially treating one material 
property as an independent variable, with other 
material properties expressed as continuous func- 
tions of the independent material property. With 
such a procedure, material selection could be treated 
as a discrete optimization variable in a nonrandom 
complex search. Similarly, the selection of com- 
mercially available elements, such as metal beams 
and columns, could be optimized by initially express- 
ing element geometric properties as continuous 
functions of a single independent variable, for example 
weight per unit length. 
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Numerous other optimization methods applicable 
to problems with nonlinear objective and constraint 
functions have been discussed in the literature, for 
example by Box [1], Hilleary [2], Schmit, Kicher, 
and Morrow [4], Schmit and. Mailed [5], and Wilde 
and Beightler [6]. The limited computational ex- 
perience reported in the presenj paper is insufficient 
to warrant a definitive claim of superiority of the non- 
random complex method over other available methods. 



The subroutines for solving the structural equa- 
tions of the elastic ring problem were coded in FOR- 
TRAN V by R. M. Slesser while he was an employee 
of the National Bureau of Standards Center for 
Computer Sciences and Technology. 



8. References 



7. Conclusion 

A nonrandom complex method has been devel- 
oped for application to optimization problems charac- 
terized by nonlinear objective and constraint func- 
tions involving continuous and/or discrete optimiza- 
tion variables. The nonrandom method appears to 
offer significant advantages over the Box [1] method 
which involves a pseudo-random process. The method 
was successful in finding approximate minimum 
weight solutions to: 

1. Seven elastic ring structural problems involv- 
ing from two to six independent variables each; and 

2. A rib-stiffened elastic plate structural problem 
involving four independent variables, two of which 
are discrete. 

These results, for problems of this complexity, sug- 
gest that the method is applicable to a large class 
of optimization problems. 



[1] Box, M. J., A new method of constrained optimization and a 

comparison with other methods, The Computer Journal 8, 

42-52 (Apr. 1965). 
[2] Hilleary, R. R., The tangent search method of constrained 

minimization, U.S. Naval Postgraduate School Tech. Rep./Res. 

Paper No. 59 (Mar. 1966). 
[3] Mitchell, R. A., Analysis of rc-degree elliptical elastic rings of 

nonuniform cross section, J. Res. NBS 72C (Engr. and Instr.). 

No. 2, 139-160 (1968). 
[4] Schmit, L. A., Jr., Kicher, T. P., and Morrow, W. M., Structural 

synthesis capability for integrally stiffened waffle plates, J. 

AIAAl,No. 12 (Dec. 1963). 
[5] Schmit, L. A., Jr., and Mallett, R. H., Structural synthesis and 

design parameter hierarchy, j. Struct. Div. ASCE 89, No. ST4, 

(Aug. 1963). 
[6] Wilde, D. J., and Beightler, C. S., Foundations of Optimization 

(Prentice-Hall, Inc., 1967). 



(Paper 72C4-282) 



258 



